#include <iostream>
#include<sys/time.h>
#include<arm_neon.h>
#define N 1000
#define num 1
using namespace std;
float A[N][N];


void reset(float A[N][N]) {
    for (int i = 0; i < N; i++) {
        for (int j = 0; j < N; j++) {
            A[i][j] = 0;
        }
        A[i][i] = 1.0;
        for (int j = i + 1; j < N; j++) {
            A[i][j] = rand();
        }
    }
    for (int k = 0; k < N; k++) {
        for (int i = k + 1; i < N; i++) {
            for (int j = 0; j < N; j++) {
                A[i][j] += A[k][j];
            }
        }
    }
}

void normal(float A[N][N]) {
    struct timeval time_start, time_end;
    float time = 0.0;
    for (int r = 0; r < num; r++) {
        reset(A);
        gettimeofday(&time_start, NULL);
        for (int k = 0; k < N; k++) {
            for (int j = k + 1; j < N; j++) {
                A[k][j] = A[k][j] / A[k][k];
            }
            A[k][k] = 1.0;
            for (int i = k + 1; i < N; i++) {
                for (int j = k + 1; j < N; j++) {
                    A[i][j] = A[i][j] - A[i][k] * A[k][j];
                }
                A[i][k] = 0;
            }
        }
        gettimeofday(&time_end, NULL);
        time += 1000 * (time_end.tv_sec - time_start.tv_sec) + 0.001 * (time_end.tv_usec - time_start.tv_usec);
    }
    cout << "普通高斯" << time / num << " ms" << endl;
}




void neon(float A[N][N]) {
    struct timeval time_start, time_end;
    float time = 0.0;
    for (int r = 0; r < num; r++) {
        reset(A);
        gettimeofday(&time_start, NULL);
        for (int k = 0; k < N; k++) {
            float* akk = new float[4];
            *akk = *(akk + 1) = *(akk + 2) = *(akk + 3) = A[k][k];
            float32x4_t vt = vld1q_f32(akk);
            delete[] akk;
            for (int j = k + 1; j + 3 < N; j += 4) {
                float32x4_t va = vld1q_f32(&A[k][j]);
                va = vdivq_f32(va, vt);
                vst1q_f32(&A[k][j], va);
            }
            for (int j = N - N % 4; j < N; j++) {
                A[k][j] = A[k][j] / A[k][k];
            }
            A[k][k] = 1.0;
            for (int i = k + 1; i < N; i++) {
                float* aik = new float[4];
                *aik = *(aik + 1) = *(aik + 2) = *(aik + 3) = A[i][k];
                float32x4_t vaik = vld1q_f32(aik);
                delete[] aik;
                for (int j = k + 1; j + 3 < N; j += 4) {
                    float32x4_t vaij = vld1q_f32(&A[i][j]);
                    float32x4_t vakj = vld1q_f32(&A[k][j]);
                    float32x4_t vx = vmulq_f32(vaik, vakj);
                    vaij = vsubq_f32(vaij, vx);
                    vst1q_f32(&A[i][j], vaij);
                }
                for (int j = N - N % 4; j < N; j++) {
                    A[i][j] = A[i][j] - A[i][k] * A[k][j];
                }
                A[i][k] = 0;
            }
        }
        gettimeofday(&time_end, NULL);
        time += 1000 * (time_end.tv_sec - time_start.tv_sec) + 0.001 * (time_end.tv_usec - time_start.tv_usec);
    }
    cout << "全部并行" << time / num << " ms" << endl;
}
//Neon第一处并行
void neon1(float A[N][N]) {
    struct timeval time_start, time_end;
    float time = 0.0;

    for (int r = 0; r < num; r++) {
        reset(A);
        gettimeofday(&time_start, NULL);
        for (int k = 0; k < N; k++) {
            float* akk = new float[4];
            *akk = *(akk + 1) = *(akk + 2) = *(akk + 3) = A[k][k];
            float32x4_t vt = vld1q_f32(akk);
            delete[] akk;
            for (int j = k + 1; j + 3 < N; j += 4) {
                float32x4_t va = vld1q_f32(&A[k][j]);
                va = vdivq_f32(va, vt);
                vst1q_f32(&A[k][j], va);
            }
            for (int j = N - N % 4; j < N; j++) {
                A[k][j] = A[k][j] / A[k][k];
            }
            A[k][k] = 1.0;
            for (int i = k + 1; i < N; i++) {
                for (int j = k + 1; j < N; j++) {
                    A[i][j] = A[i][j] - A[i][k] * A[k][j];
                }
                A[i][k] = 0;
            }
        }
        gettimeofday(&time_end, NULL);
        time += 1000 * (time_end.tv_sec - time_start.tv_sec) + 0.001 * (time_end.tv_usec - time_start.tv_usec);
    }
    cout << "第一处并行 " << time / num << " ms" << endl;
}
//Neon第二处并行
void neon2(float A[N][N]) {
    struct timeval time_start, time_end;
    float time = 0.0;
    for (int r = 0; r < num; r++) {
        reset(A);
        gettimeofday(&time_start, NULL);
        for (int k = 0; k < N; k++) {
            for (int j = k + 1; j < N; j++) {
                A[k][j] = A[k][j] / A[k][k];
            }
            A[k][k] = 1.0;
            for (int i = k + 1; i < N; i++) {
                float* aik = new float[4];
                *aik = *(aik + 1) = *(aik + 2) = *(aik + 3) = A[i][k];
                float32x4_t vaik = vld1q_f32(aik);
                delete[] aik;
                for (int j = k + 1; j + 3 < N; j += 4) {
                    float32x4_t vaij = vld1q_f32(&A[i][j]);
                    float32x4_t vakj = vld1q_f32(&A[k][j]);
                    float32x4_t vx = vmulq_f32(vaik, vakj);
                    vaij = vsubq_f32(vaij, vx);
                    vst1q_f32(&A[i][j], vaij);
                }
                for (int j = N - N % 4; j < N; j++) {
                    A[i][j] = A[i][j] - A[i][k] * A[k][j];
                }
                A[i][k] = 0;
            }
        }
        gettimeofday(&time_end, NULL);
        time += 1000 * (time_end.tv_sec - time_start.tv_sec) + 0.001 * (time_end.tv_usec - time_start.tv_usec);
    }
    cout << "第二处并行 " << time / num << " ms" << endl;
}


int main()
{

    normal(A);

    neon(A);

    neon1(A);

    neon2(A);


    cout << endl;
    return 0;
}